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In this grant period, the focus has been on enhancement and 
application of the direct simulation Monte Carlo (DSMC) particle 
method for computing hypersonic flows of re-entry vehicles. 
Enhancement efforts dealt with modeling gas-gas interactions for 
thermal non-equilibrium relaxation processes and gas-surface 
interactions for prediction of vehicle surface temperatures. Both 
are important for application to problems of engineering interest. 
The code was employed in a parametric study to improve future 
applications , and in simulations of aeropass maneuvers in support 
of the Magellan mission. 

Detailed comparisons between continuum models for internal energy 
relaxation and DSMC models reveals that several discrepancies 
exist . These include definitions of relaxation parameters and the 
methodologies for implementing them in DSMC codes. These issues 
were clarified and all differences were rectified in a paper 
(Appendix A) submitted to Physics of Fluids A, featuring several 
key figures in the DSMC community as co-authors and B. Haas as 
first author. This material will be presented at the Fluid 
Dynamics meeting of the American Physical Society on November 21, 
1993. 


The aerodynamics of space vehicles in highly rarefied flows are 
very sensitive to the vehicle surface temperatures. Rather than 
require prescribed temperature estimates for spacecraft as is 
typically done in DSMC methods, a new technique was developed 
which couples the dynamic surface heat transfer characteristics 
into the DSMC flow simulation code to compute surface temperatures 
directly. This model, when applied to thin planar bodies such as 
solar panels, was described in AIAA Paper No. 93-2765 (Appendix B) 
and was presented at the Thermophysics Conference in July 1993. 



The paper has been submitted to the Journal of The rmophy sics and 
Heat Transfer . 

Application of the DSMC method to problems of practical interest 
requires a trade off between solution accuracy and computational 
expense and limitations. A parametric study was performed and 
reported in AIAA Paper No. 93-2806 (Appendix C) which assessed the 
accuracy penalties associated with simulations of varying grid 
resolution and flow domain size. The paper was also presented at 
the Thermophysics Conference and will be submitted to the journal 
shortly . 

Finally, the DSMC code was employed to assess the pitch, yaw, and 
roll aerodynamics of the Magellan spacecraft during entry into the 
Venus atmosphere at off-design attitudes. This work was in 
support of the Magellan aerobraking maneuver of May 25-Aug. 3, 
1993. Furthermore, analysis of the roll characteristics of the 
configuration with canted solar panels was performed in support of 
the proposed "Windmill" experiment. Results were reported in AIAA 
Paper No. 93-3676 (Appendix D) presented at the Atmospheric Flight 
Mechanics Conference in August 1993, and were submitted to Journal 
of Spacecraft and Rockets. 



Appendix A 
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ABSTRACT 

For internal energy relaxation in rarefied gas mixtures, exact relationships are derived between the 
selection probability P employed in direct simulation Monte Carlo (DSMC) methods and the macroscopic 
relaxation behavior dictated from Jeans’ equation through collision number Z . These expressions apply 
to the Borgnakke- Larsen model for internal energy exchange mechanics. Simulation results prove that the 
common assumption P = 1 /Z, although it can be a fair approximation, is invalid in general since relaxation 
rates depend upon the intermolecular potential, the number of internal degrees of freedom, and the selection 
methodology. Differing definitions of collision number Z appearing in the literature are clarified. As proven 
theoretically, all relaxation curves coalesce to a single curve when plotted against the cumulative number of 
collisions, regardless of the intermolecular potential. These same relationships apply for variable Z as well. 
The present work concludes that DSMC kinetics in gas mixtures only match Jeans- Landau-Teller behavior 
exactly when prohibiting multiple relaxation events during a single collision. 
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I. INTRODUCTION 

Direct simulation Monte Carlo (DSMC) 1 ’ 2 and its variants 3-8 have proven to be reliable tools in 
computing rarefied flows. These methods employ thousands or millions of discrete computational particles 
whose motion and interaction simulate gas dynamics. High temperature, rarefied flows are often char- 
acterized by thermal nonequilibrium in which various molecular energy modes do not share a common 
Boltzmann distribution of energies. Simulating these flows, and the corresponding relaxation gas dynamics, 
requires application of simple phenomenological models which capture the essential macroscopic behavior 
while retaining computational efficiency in DSMC methods. Several models have appeared in the litera- 
ture addressing rotational, 9- 11 vibrational, 12-15 and chemical 13 ’ 16-20 nonequilibrium. A common element 
within most of these models is application of the Borgnakke-Larsen method 21 in distributing post-collision 
thermal energy among the energy modes of a pair of colliding particles. This technique assigns internal 
and translational energies by sampling from the corresponding equilibrium distributions. Of vital concern, 
however, is the resulting macroscopic rate at which this technique promotes relaxation of the ensemble gas 
toward equilibrium. 

DSMC methods test colliding particles for possible internal energy exchanges or chemical reactions. 
The collision selection rules employ probabilities which must be defined so as to result in desired macroscopic 
relaxation rate behavior. This typically involves identifying a macroscopic rate equation which the simulation 
attempts to match through appropriate definitions of selection probabilities P. These rate equations often 
contain parameters such as collision numbers Z , characteristic times r, or rate coefficients which are estimated 
from experiments or analysis. Ambiguity regarding specific definitions of the targeted macroscopic rates and 
coefficients has led to confusion in defining appropriate selection probabilities in DSMC methods. This is 
particularly evident in studies of rotational relaxation for which the Borgnakke-Larsen collision mechanics 
are widely used and the relaxation rates may be clearly derived. 

In studies of vibrational relaxation, Haas 13 developed a relationship between the appropriate DSMC 
selection probability P and the resulting Landau-Teller macroscopic relaxation behavior. This relationship 
was later generalized, allowing for rotational relaxation, by Lumpkin et a/. 22 in the form of a “correction 
factor” for collision numbers Z. It was observed that the DSMC probabilities P in general cannot be specified 
by the assumption P = 1/Z y which has been employed throughout the DSMC community. Dissimilar 
applications of the selection probability and incompatible definitions of collision number led to disagreement 
among DSMC researchers regarding these correction factors. 

The present work clarifies the definitions of macroscopic rate parameters employed in rotational re- 
laxation, including differing interpretations of Jeans’ equation. 23 It will be shown that relaxing energies, 
when plotted against the cumulative number of molecular collisions, lead to a single solution regardless of 
the assumed intermolecular potential. Differing definitions of the rotational collision number Z and colli- 
sion time t c will be clarified. The two most widely used selection methodologies are discussed, and DSMC 
rotational relaxation probabilities P specific to each methodology will be derived to match the macroscopic 
rate equations exactly. Prevalent selection methods cannot reproduce Jeans’ equation exactly for the case 
of heteromolecular collisions, nor can they reproduce Jeans and Laundau-Teller 24 rate behavior for simul- 
taneous rotational and vibrational relaxation. Thus, a new selection methodology is introduced that is 
better suited to the general case of gas mixtures experiencing multimode internal energy relaxation. Finally, 
considerations regarding variable Z will be addressed. 
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II. DEFINITIONS OF RELAXATION TIMES AND COLLISION NUMBERS 

Most continuum analyses and computational codes employ some form of the Jeans- Laundau-Teller 
equation 23,24 for study of internal energy relaxation in gases. For rotational relaxation in an single-species 
adiabatic reservoir, this equation employs the instantaneous equilibrium rotational energy E*(t) and a time 
constant given by the product of the instantaneous collision time r c and a characteristic collision number Z 
as follows: 


dEr _ E*{t) ~ Erjt) 

dt T e Z W 

Here, E r is the mean rotational energy per molecule in the ensemble, and t represents time. Note that E?(t) 
varies proportionally with translational temperature T, which will itself vary in time during relaxation. 
Similarly, r c may vary with T, depending upon the intermolecular potential assumed for the gas. Most 
DSMC methods employ the Variable Hard Sphere (VHS) model of Bird 25 which assumes rigid-sphere col- 
lision mechanics but establishes the collision rate in a manner consistent with an inverse-power potential 
with exponent a. This parameter varies between the limits for Maxwell molecules (a = 4) and hard sphere 
molecules (a = oc). Employing a collision cross-section of the form where g is the relative trans- 

lational speed of the pair and c is a constant, DSMC methods lead to an exact definition of collision time 
given by 


T c ~ 



( 2 ) 


Here m* is the reduced mass of the pair, n is the number density, and k is Boltzmann’s constant. Note 
that r c is independent of T for Maxwell molecules. 

An alternative expression for macroscopic relaxation may be expressed in terms of the final equilibrium 
rotational energy E*( oo) and a new collision number Z as follows: 


dEr. = E*(oo)-Er(t) f3) 

dt r c Z 

In order to relate Z and Z, it is noted that at any time t for a single species adiabatic reservoir, the first law 
of thermodynamics yields 


Et(t) + E r (t) = E 0 (4) 

where E t is the mean translational energy and E 0 is the mean total energy which is constant. Moreover, 
the instantaneous equilibrium mean rotational energy is related to E t by definition of the thermal DOF for 
rotation £ r (=2 for diatomic species) and translation £* = 3 as follows, 

£*(<) = ^E t {i). (5) 

By virtue of Eqs. (4) and (5), the final equilibrium rotational energy is expressed by 

Cr 

Cl + Cr 


E?(oo) = 


E, 


( 6 ) 
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E:(t)-E r (t)= (l + l) 
Comparing Eqs. (1), (3), and (7), one finds 


£;*(oo) - E r (t ) 


( 7 ) 


l = ioTz) z - (8) 

Both relaxation equations may be integrated exactly for Maxwell molecules, assuming contant collision 
numbers and densities in the adiabatic reservoir, leading to the following solutions 


1 - 


Erjt) 

E*(oo) 


(~iz nr 1 ) = (~zi) ■ 


= exp ( l = exp [ J . (9) 

For all other types of molecular interaction, where r c is dependent upon temperature, numerical integration 
is required. Note that the relation expressed in Eq. (8) holds irrespective of the interaction potential; i.e., it 
is not limited to the case of Maxwell molecules. 

As shown in Fig. la for particular values of a, both rate equations lead to the same relaxation behavior 
in the adiabatic reservoir when the collision numbers are related as suggested in Eq. (9). Furthermore, 
these solutions coalesce to a single relaxation curve in Fig. lb, regardless of VHS parameter a, when plotted 
against a time scale t representing the cumulative number of collisions per molecule given analytically by 
the expression, 


t = 



dt_ 

r c 


( 10 ) 


Note that the curve is purely exponential in t ) relaxing to 1/e at t/Z = 1. 

The conclusion here is that Eqs. (1) and (3) reflect the same relaxation behavior but simply employ 
differing definitions of the collision number as noted in Eq. (8). 


A. Definitions of mean collision time 

Since collision numbers are often estimated from observed or computed relaxation behavior, one must 
ensure that consistent definitions of mean collision time r c are employed to reduce the Z data. 26 Many 
experimental definitions 27 of collision time r c are quoted assuming either the Chapman-Enskog (C-E) first 
approximation for viscosity, or the analogous result derived from simple kinetic theory (KT), leading to the 
respective expressions, 


■ _ 16 / 5 /» | _ TW 

'c[ C _E — t 1 C j KT 4P’ 


( 11 ) 


where P = nkT is pressure and fi is viscosity estimated from experiments. Recall, however, that r c is defined 
by the ratio A/(c), where A is the mean free path and (c) = (SkT/nm) 1 / 2 is the mean thermal speed for 
molecular mass m. Since A has unambiguous meaning for VHS molecules regardless of a used in the DSMC 
technique, Bird 28 determined a clear relationship between A and /i by employing the Chapman-Enskog first 
approximation, leading to: 


e ' VHS 30 P‘ 


( 12 ) 
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general. As a result of the quadratic term in Eq. (14), there exists a minimum allowable value of Z to insure 
that P PART is a real number less than or equal to unity, 

f 2, for diatomic hard spheres; 

> \ l /3, for diatomic Maxwell molecules. 

To evaluate the models presented here, adiabatic reservoirs composed of a single diatomic species were 
simulated during rotational relaxation for VHS parameters a = {4, oo}. Each reservoir was initialized with all 
energy in the translational mode only; i.e., E r { 0) = 0. Both selection methodologies were used with a targeted 
relaxation rate corresponding to Z - 10. As plotted in Figs. 3a and 4a, results for simulations employing the 
assumption P —\jZ differ from the theoretical curve of Eq. (9) and from one another for particle-selection 
and pair-selection methodologies, respectively. These differences not only prove the assumption is faulty, 
but they also demonstrate the dependence of the relaxation rate upon VHS parameter a. Repeating these 
simulations, using the probability computed from Eqs. (14) and (13), led to the plots in Figs. 3b and 4b 
which all collapse upon the theoretical curve correctly. The two figures appearing in Ref. (22) demonstrated 
similar results. Unfortunately, due to a typographical error, the figure legend and caption from each graph 
were switched with those of the other. 

While the pair-selection method enjoys simplicity and computational efficiency, it is not well suited 
to general gas mixtures involving differing relaxation rates and internal DOF per species. These same 
difficulties apply for the particle-selection method unless one prohibits the possiblility of both particles 
relaxing internally during a single binary collision (double relaxation). To explain this problem, note that 
for rotational relaxation of species A due to collisions with molecules of species B, Jeans' equation is expressed 
as 


dE rA 

dt 


E; A (t)-E rA (t) ' 

t Ca\b Za\B 


(15) 


Notice that this expression is independent of the rotational energy E re of the partner species. However, 
for the pair-selection method in the case of heteromolecular collisions, relaxation of A molecules would be 
influenced by the rotational energy of B molecules since both relax simultaneously. This is also true for the 
particle-selection method because a possibility exists that both particles will relax during one collision. As a 
consequence, E Tb would appear in the resulting relaxation rate equations for species A such that it cannot 
be cast in the form given in Eq. (15) for heteromolecular collisions. 

A similar problem arises when coupling rotational and vibrational relaxation. If collisions permit 
simultaneous relaxation of each mode, the corresponding rate equations for each cannot be cast into the 
Jeans-Landau-Teller form exactly because vibrational energies would appear in the rotational relaxation 
rate equation and vice-versa. As derived at the end of the Appendix, however, decoupling the relaxation of 
each mode does lead to relaxation behavior dictated in the Jeans-Landau-Teller form. 

Nonetheless, the methods above work very well for vibrationally-frozen gases in which the rotational 
energies of colliding species have the same number of DOF and relax at the same rate. Such assumptions are 
employed frequently in continuum and DSMC computations for real gas flows, although the present work is 
not limited to these assumptions. 
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B. General DSMC selection method 

The intent of the following is to introduce a new selection methodology which is not restricted to the 
assumptions above and is therefore applicable to gas mixtures. The method is adapted from the particle- 
selection method in that probabilities P are applied to each molecule individually within the colliding pair, 
but the event in which both particles relax is prohibited (called particle-selection without double relaxation). 
This insures that relaxation of A molecules is independent of the energy of B molecules, and the Borgnakke- 
Larsen mechanics can therefore reproduce the relaxation rate behavior defined in Eq. (15). The appropriate 
functional form for the relaxation probability P f pART is derived in the Appendix for particle-selection prohibit- 
ing double-relaxation in the case of binary collisions in general gas mixtures. When applied to a single-species 
gas, this relationship is given by 



where Z was defined in Eq. (1) and is related to Z in Eq. (8). Note that P pART depends upon the VHS 
parameter as well as the number of internal DOF experiencing relaxation. Again, restrictions on allowable 
values Z are imposed by the quadratic term; 

2 f 3, for diatomic hard spheres; 

^ ( 10/3, for diatomic Maxwell molecules. 

As before, the method is evaluated through the use of an adiabatic reservoir of a single-species diatomic 
gas. Figure 5 shows relaxation behavior with application of uncorrected and corrected probabilities with the 
targeted relaxation collision number Z = 10. Note that the assumption P = \/Z is close but not exact. 
Nonetheless, specific relationships between Z and P may be derived to yield the correct relaxation rate 
behavior in a single species gas regardless of the selection methodology employed in the DSMC technique. 
The appropriate function, however, is unique for each particular selection methodology as demonstrated by 
the differences between Eqs. (13), (14), and (16). Depending upon the methodology employed, however, the 
approximation P = 1 f Z may be fairly accurate under practical flow conditions. 

C. Variable relaxation probabilities 

Although the results noted above correspond to invariant relaxation rates Z, the expressions for P 
are not restricted to that assumption. Parker 30 expressed the dependence of Z upon temperature T in a 
homonuclear diatomic gas. Similar observations were made analytically by Lordi and Mates 31 and exper- 
imentally by Carnevale et al 27 For DSMC methods, Parker’s expression may be used to evaluate a single 
value of Z for all collisions within a computational cell of local translational temperature T. It can, however, 
be computationally expensive to calculate T which may itself be ill-defined in nonequilibrium flows. 

Alternatively, Boyd 9 developed an expression for Z which is dependent upon the relative translational 
energy of individual collisions and which reproduces Parker’s expression when integrated over all collisions at 
equilibrium. However, application of probabilities computed in this manner lead to energy distributions that 
are biased away from the Boltzmann form assumed in the Borgnakke-Larsen model for relaxation mechanics. 
As such, the model will fail to achieve detailed balance. Boyd countered this by applying a single “averaged” 
probability throughout the cell. Abe 10 achieved detailed balance while retaining individual collision-specific 
probabilities by incorporating the biased distributions in his adaption of the Borgnakke-Larsen model for 
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relaxation mechanics. Alternatively, one could employ Z as a function of the sum of translational and 
rotational energies since this would not bias the distributions normally employed in the Borgnakke-Larsen 
technique to partition the energy into individual modes. The important point is that regardless of the 
formulation used to determine Z, its relationship to P still applies as defined in the present work. 

V. CONCLUDING REMARKS 

Distinct relationships exist between the microscopic probabilities P employed at the collision level in 
DSMC methods and the resulting macroscopic relaxation rate behavior in the gas as described by collision 
numbers Z. These relationships are sensitive to the selection methodology; i.e., the particular manner in 
which the probabilities are applied. For relaxation of internal modes through application of the Borgnakke- 
Larsen model of energy-exchange mechanics, the selection probability is not given generally by the simple 
assumption P = 1/Z. Indeed, it is dependent upon the intermolecular VHS potential and the number of 
DOF contributing to the relaxation process. Furthermore, there are differing definitions appearing in the 
literature regarding collision numbers Z and collision times r c . In most cases a simple relationship resolves 
these differing definitions as identified in the present work. It was also observed that regardless of the VHS 
potential employed, all theoretical relaxation curves coalesce to a single curve when plotted against a time axis 
representing the cumulative number of molecular collisions. Simulations proved that the DSMC technique 
leads to the same relaxation behavior as dictated by the macroscopic rate equation only when employing 
the probabilities derived in the present work, although P = 1/Z can be a fair approximation in practical 
scenarios when employing the particle-selection methodology. Finally, DSMC kinetics in gas mixtures only 
match Jeans-Landau-Teller behavior exactly when prohibiting multiple relaxation events during a single 
collision. 
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APPENDIX: DERIVATION OF RELAXATION PROBABILITY 

The relationship between DSMC probabilities and macroscopic collision numbers can be derived by 
assessing the effect of relaxation collisions upon the internal energy of the entire ensemble of particles, 
as a result of Borgnakke-Larsen mechanics, and casting the resulting expression into the form dictated in 
Eq. (1). The following development applies to the particle-selection methodology without double-relaxation, 
a technique most readily extended for general gas mixtures. The derivation for particle-selection including 
double-relaxation is similar but may only be applied to gases for which all species enjoy identical relaxation 
behavior. Under the same limitations, the derivation for pair-selection appears in Lumpkin et al. 

A • Fraction of relaxing particles 

Consider a gas mixture for which the macroscopic rotational relaxation of species A due to collisions 
with particles of species B is governed by collision number Z\ \b in Jeans’ Eq. (15). Likewise, relaxation of 
species B due to collisions with species A is governed by Z b \a • In this mixture, the total relaxation rate of 
species A is given by the sum of contributions resulting from collisions with all possible partner species B in 
the gas, 




(Al) 



E r *{t)-E rA {t) . ^ d±r >t 
T 'a |B z a\b 



In many continuum analyses, the sum over binary collision time constants t Ca ^ b Z a \b is often replaced by 
a molar average over mean collision times; Eq. (Al) is employed in the present work as it is simpler and 
more general. Note that a new notation is employed in this appendix to facilitate derivation for general gas 
mixtures. We wish to derive corresponding expressions for DSMC probabilities P a \b and Pb\a- 

In a given collision, the first particle of the pair is tested for relaxation. If accepted, the Borgnakke- 
Larsen method for collision mechanics is applied for relaxation of this particle only. If rejected, then and 
only then will the second particle be tested for relaxation which may or may not be accepted. This selection 
methodology is depicted in the probability tree of Fig. 6. Note that the probability of the first particle being 
of species A in a given A+B pair is 1/2, assuming random pairing. From the figure, it is clear that the fraction 
of collisions for which particle A is the first particle in the pair, and is accepted for relaxation, is Pa\b/2> The 
fraction of collisions for which A is the second particle, yet is accepted for relaxaton, is (1 — Pb\a)Pa\b /%- 
Note also that the fraction of accepted collisions equals the fraction of accepted particles for this selection 
methodology. It follows that the total fraction F a \b of A particles which collide with B particles and 
experience relaxation is given by the sum 


( 


Fa\b = \Pa\b + \ (l - p b\a) Pa 


— Pa\b 



(A2) 

(A3) 


B. Effect upon ensemble energy 

Denote the mean rotational energy of the ensemble of A particles at a given time t in the DSMC 
simulation by E* s E rA (t). In a relaxing A + B collision, the mechanics of Borgnakke- Larsen will lead to 
an equilibrium division of total mean collision energy E x t 4* where E X T is the mean relative translational 
energy of colliding particles at time t and has C T = 4 — 4/a DOF near equilibrium. In this model the 
post- collision rotational energy of those particles at time t which relax during the time step (denoted E* a 
is given by the ratio of equilibrium DOF, 22 


El' = k* . (E‘ + E* ) , (A4) 

r A,tlx f _J_ £ V T r A / * V ' 

v> 

where is the equilibrium number of rotational DOF for species A. 

During a given time step of duration At, the total fraction of species A which collide with species B is 


given by the ratio At/r CA ^ B . Here, t Ca | b is the collision time for species A given in Eq. (2) when using n = n B 
and a = a^ B . The mean rotational energy at the end of the time step (denoted would include 

contributions from relaxing particles and from particles unchanged during the time step, 
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Substituting the expression from Eq. (A4) and simplifying leads to 


= 


At 


' c a\ b 


Fa\b 


<T 

Cr + C r A 


V 



(A6) 


C. Relaxation rate equation 

To cast the result of Eq. (A6) into a form representative of Jeans Eq. (15), we note that the DSMC 
method approximates the rate of rotational energy change as, 


jga - - Sa 

dt ~ At 

Cr UrMK-Pr* 

<T + Cr„ 

If the ensemble gas is near equilibrium in the translational mode such that molecular velocities are 
distributed in Maxwellian form, then the relationship between mean translational energy E t and mean 
relative translational energy of colliding particles E r is given by the ratio of DOF , 


(A7) 

(A8) 


EL 

E ‘r 


Cr' 


(A9) 


Likewise, the instantaneous equilibrium rotational energy E* (t) is defined by the expression, 


E: A (t) = ^ E l (A10) 

Substituting Eqs. (A9) and (A 10) into Eq. (A8) leads to the rate equation 

jgL _ ' Cr 

dt C r + Cr A 

^By-comparisomtorBqr-(15), and employing the expression in Eq. (A2), the following relationship results 
between the macroscopic collision number and the required DSMC probabilities, 






( t 'a\b/ F *\b) 


(All) 


1 

%a\b 


= Pa\b 



Cr 
C T + 


(A12) 


Note that Z a \b is dependent upon both P A \ B and P B \ A ■ Substitution of the appropriate indices leads to 
the corresponding expression for Z B \ A . When combined, these equations yield a single quadratic expression 
which may be solved for P A \ B in terms of parameters associated with both species A and B, 


lp2 _ ( Kr+Cr, 1 

2 \2 C T Z A \ B 


1 Cr+_C^._L_ + 

2 Cr Z B \ A ) 


Pa\b + 


C T + Cr A 1 
C T Z A \ B 


= 0. 


(A13) 
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A similar relation for Pb\a is formed by switching the appropriate indices. The special case of a single- 
species gas, corresponding to Q rA — Cr B and Za\b = %b\a> leads to the result cited in Eq. (16). Similarly, 
collision between a polyatomic molecule (£ r>t > 2) and an atom (£r s = 0, Zb\a ~ °°) Is equivalent to the 
pair-selection result cited in Eq. (13). 

To demonstrate application of this selection method for a general binary mixture, relaxation of an 
adiabatic reservoir of species A and B was simulated where all energy was initially in the translational 
mode only. The reservoir was specified a s follows: n A /n = 0.3, n B /n = 0.7, a = co, ( \ rA — 2, £r B = 3. 
The relaxation rates and corresponding values of Z, Z, and P are specified in Table 1. Results of DSMC 
simulation permitting double relaxation and using P = 1/Z are plotted in Fig. 7 along with results when 
prohibiting double-relaxation and using P defined in Eq. (A13). Note that the latter results compare best to 
those from numerical integration of Jeans Eq. (Al). This serves to validate the selection rule of the present 
work for gas mixtures. 

The only limitation to the analysis above is the assumption in Eq. (A9) that the translational mode is 
near equilibrium. Rotational energy exchanges which are too frequent during one time step in a relaxing gas 
may lead to an insufficient number of intermediate elastic collisions, thus preventing the translational mode 
from equilibrating locally. As a consequence, the assumption in Eq. (A9) may be invalid when the rotational 
selection probability Pa\b approaches unity. 

This derivation applies only to the particle-selection methodology which prohibits double-relaxation 
on one collision. For the methodology which permits double relaxation, or for the pair-selection method, 
the corresponding expressions for could not be reduced to the simple form in Eq. (All) as required 

to match Jeans’ equation for general gas mixtures. For simple gases in which the rotational modes of all 
species are assumed equal, the derivation above may be tailored to these other methodologies, leading to 
the expressions in Eqs. (13) and (14). In those cases, however, one must be careful to employ the correct 
particle-fractions which may differ from the corresponding collision-fractions. 


D. Multimode internal relaxation 

If collisions permit simultaneous relaxation of different internal energy modes of a molecule, then the 
energy relaxation rate equation corresponding to Borgnakke-Larsen mechanics for one mode would include 
energy terms for the other mode. As such the rate equations could not be cast into the Jeans-Landau-Teller 
form in Eq. (15). The macroscopic behavior of coupled multimode relaxation is not well understood, and 
perhaps may not be described correctly by Eq. (15). Nonetheless, to reproduce this macroscopic behavior 
exactly, one needs simply to prohibit simultaneous relaxation of two energy modes during a single collision. 
One such selection methodology is depicted in a probability tree in Fig. 8. 

Because rotational relaxation of each particle is tested prior to testing for vibrational relaxation in 
Fig. 8, the resulting rotational probabilities are the same as derived above and described in Eq. (A13). 
However, expressions for vibrational relaxation probabilities P v include rotational probability terms P r and 
were derived in the same manner used above, leading to 
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Here, t„ a ^ b refers to the characteristic vibrational relaxation time in the Landau-Teller rate equation and 
has been measured experimentally by Millikan and White 32 and Camac 33 for several gases. The effective 
DOF <„ reflects the vibrational energy content and is dependent upon flowfield temperature. 14 - 24 

Switching the appropriate indices in Eq. (A14) leads to the corresponding expression for P Vb \ a - Note 
that the selection order in the methodology of Fig. 8 could be rearranged freely, leading to different probability 
expressions which are nonetheless similar to those in Eqs. (A13) and (A14). 
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TABLE 1: Relaxation rate parameters for simulated mixture 


i\j a Z Z 1/P b 


A l 

A 

6 

3.6 

3.4142 

Ai 

B 

8 

4.8 

4.9876 

B j 

A 

15 

7.5 

7.7122 

B 

B 

20 

10.0 

10.9045 


Denotes relaxation of species i via collisions with species j. 
^ From Eq. (A13). 
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FIGURE CAPTIONS 

Fig. 1 Solutions of Jeans’ Eqs (1) and (3) for hard-sphere and maxwell molecules, using collision numbers 
related by Eq. (8). a) Results plotted against time normalized by the initial t c . b) Results plotted against 
the cumulative number of collisions t defined in Eq. (10). 

Fig. 2 Ratio of P from Eq. (14) to 1/Z, derived for several VHS intermolecular potentials, a, corresponding 
to the particle-selection methodologies for which double relaxation is permitted in a diatomic gas « r = 2). 
Fig. 3 Diatomic rotational energy during relaxation in an adiabatic reservoir ( Z = 10). DSMC results 
used the particle-selection methodology permitting double-relaxation with P defined from: a) 1 /Z, and 
b) Eq. (14). 

Fig. 4 Diatomic rotational energy during relaxation in an adiabatic reservoir ( Z = 10). DSMC results used 
the pair-selection methodology with P defined from: a) 1 /Z, and b) Eq. (13). 

Fig. 5 Diatomic rotational energy during relaxation in an adiabatic reservoir (Z = 10). DSMC results 
used the particle-selection methodology prohibiting double-relaxation with P defined from: a) 1 /Z t and 
b) Eq. (16). 

Fig. 6 Probability decision tree for the particle-selection methodology prohibiting double-relaxation. Joint 
probabilities for all possible outcomes are determined on the far right. 

Fig. 7 Rotational energy during relaxation of an adiabatic reservoir of the binary gas mixture described 
in Table 1; comparison between solution of Jeans Eq. (Al) and DSMC simulation, a) permitting double 
relaxation and using P = 1 /Z, and b) prohibiting double relaxation and using P defined from Eq. (A13). 
The time axis is normalized by the initial value of t Ca | a . 

Fig. 8 Probability decision tree for rotational and vibrational relaxation using a particle-selection method- 
ology which prohibits multiple relaxation events per collision. 
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Models for Dynamic Surface Temperatures 
During Rarefied Aeropass Maneuvers 

Brian L. Haas* 

Eloret Institute , Palo Alto , California 94303 

Models are developed to compute surface temperatures in a direct simulation Monte Carlo ( DSMC ) 
particle method for vehicles during atmospheric entry. This involves time-integration of the governing 
heat transport equations over the duration of the time-dependent heating pulse , and then coupling these 
temperatures into the DSMC solution. This model is adapted specifically for application to thin panels in 
highly-rarefied flows and yields the transient thermal response of individual surface facets by accounting 
for radiative emission, thermal conducti vity and material heat capacity The model may be simpified for 
materials with instant thermal response due to negligible heat capacity Simulation results for radiating 
panels with internal conductivity are validated by comparison to corresponding free molecule theory ; 
These models are robust , contribute negligible additional computational burden , and permit calculation 
of all surface-facet temperatures independently. Furthermore , these models may be used in continuum 
methods for computational fluid dynamics. Sensitivity of aerodynamic heating , drag , and lift to surface 
temperature is assessed through simulations of several hypersonic rarefied flows past hat panels at 45° 
incidence . 


Nomenclature 

A thermal accommodation coefficient 

C aerodynamic force coefficients 

c effective material heat capacity (J/m 2 K) 

E error residuals for Newton iteration 

H flow similarity parameter, Eq. (14) 

k Boltzmann constant, 1.3805 x 10 -23 J/K 

Kn Knudsen Number 

M Mach number 

n gas number-density (m~ 3 ) 

q transient net convective heat flux 

Q net convective heat flux computed in DSMC 
S speed ratio, S = \Zi/2M 00 

T temperature 

T temperature normalization parameter, Eq. (4) 
t time during heating pulse 

Z number flux parameter, Eq. (17) 

a exponent of intermolecular potential 

P angle of incidence of panel to free-stream flow 
7 ratio of specific heats 

e material radiative emissivity 

C int molecular internal degrees of freedom (£ im% = 2) 
i) temperature similarity parameter, Eq. (15) 

k effective thermal conductivity (W/m 2 K) 

p gas mass-density (kg/m 3 ) 

<t Stefan-Boltzmann const., 5.67 x 10~ 8 W/(m 2 K 4 ) 
r heating pulse time-constant, perigee at t = t 

u molecule translational energy-flux, Eq. (16) 

* Research Scientist, Aerothermodynamics Branch. Member, AIAA. 
Mailing Address: NASA Ames Research Center, MTS 230-2, 
Moffett Field, CA 94035- 1 000. 

This paper is declared a work of the UJ5. Go vernm ent and is not 
subject to copyright protection in the United Stale*. 


Subscripts: 

a , b front,back-side of panel 

d drag or deep-space value 

l lift 

o stagnation value 

i, j value at surface node, neighboring node 

w value at wall or surface 

o© freestream value 

Introduction 

To perform orbital maneuvers in a manner which ex- 
pends little onboard fuel, spacecraft may pass through the 
outer reaches of planetary atmospheres. Such scenarios are 
characterized by extremely rarefied flows at high velocities 
with transient heat fluxes to vehicle surfaces. However, del- 
icate spacecraft components such as solar panels may not 
withstand severe aerodynamic heating or body loading. Ac- 
curate prediction of the rarefied flow field and subsequent 
vehicle aerodynamics is therefore essential to successful 
mission planning for aeropass maneuvers. One contempo- 
rary example of such a scenario is the aerobraking of the 
Magellan spacecraft to circularize its highly eccentric orbit 
about Venus. 1 Magellan’s solar panels, which are mounted 
normal to the flow direction, place severe restrictions on 
the acceptable flight altitudes and body torques during each 
aeropass. More importantly, non-zero material heat capac- 
ity and finite conductivity through the panels lead to tran- 
sient surface temperatures during the short heat pulse of aer- 
obraking. 

Highly rarefied flows, for which the ratio of molecular 
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mean free path to a body dimension is large (Knudsen num- 
ber, Kn>0.10), are best simulated computationally with di- 
rect simulation Monte Carlo (DSMC) particle methods. 2 3 
Here, gas dynamics is modeled directly by the motion and 
interaction of thousands or millions of discrete particles. 
Particles which strike the vehicle may reflect back into the 
flow with velocities and internal enemies corresponding to 
full or partial accommodation to surface. In rarefied flows, 
molecular diffusion from vehicle surfaces affects the in- 
coming flow. Therefore, surface temperatures play a signif- 
icant role in establishing the net heat transfer to the body. 
More importantly, aerodynamic body forces perpendicular 
to the freestream direction depend only upon the flux of 
particles away from the surface. Since this efflux is itself 
greatly dependent upon surface temperature, the net aero- 
dynamics of a vehicle in rarefied flows exhibits tremendous 
sensitivity to surface temperatures as well. 

The objectives of the present work are two-fold. The 
first goal is to assess the sensitivities of aerodynamic heat- 
ing and forces to surface temperature for several blunt-body 
entry flows. The second goal is to develop new models 
which calculate dynamic surface temperatures directly in 
the particle simulation rather than requiring the use of pre- 
scribed isothermal estimates. These models account for the 
radiative emissivity, thermal conductivity, and heat capac- 
ity of the material. These heat-transport properties are par- 
ticularly important for atmospheric entry heating pulses of 
short duration in which transient temperature profiles lead 
to peak temperatures occurring after the peak heating point 
because of non-zero heat capacity of the material. Simplif- 
ing the dynamic model for materials with negligible heat 
capacity permits validation of the methods through com- 
parison to theoretical predictions in free molecular flows. 
Although models of this type have been demonstrated for 
application to ordinary single-sided surfaces, 4 the present 
study will address application to thin panels which have dy- 
namic temperatures on both front and back sides. 


= {5,10,20} 
T./T" ={1,2,4,8,16} 


In these simulations, the wall temperature was fixed at 
prescribed value T w over the entire plate surface. Reflec- 
tion from the surface was fully diffuse and thermally ac- 
commodated (A = 1 .0). Figure 2 depicts the dependence 
of the steady-state net convective heat transfer Q upon T w . 
Note that the ratio QiT^/QiT^) is independent of Kn^ . 
Note also that the curves are nearly identical at all Mach 
numbers when surface temperature is normalized by the to- 
tal (stagnation) temperature, T 0 . Dependence of heating 
upon M^ is therefore described by the standard isentropic 
expression for T 0 , 



It follows that the dependence of Q upon T w dimin- 
ishes with increased M^ , regardless of flow-field rarefac- 
tion, Kji . 

Sensitivity of plate drag to surface temperature is 
demonstrated in Fig. 3. Just as with heating, drag sensi- 
tivity is nearly identical for all Kn^ , but depends strongly 
upon M^ . In general, drag is marginally sensitive to wall 
temperature, changing by less than a factor of two over the 
range of T w studied here. 

Lift coefficients are more sensitive to T w as demon- 
strated in Fig. 4. Unlike heating and drag, lift depends 
greatly upon Kn^. For example, at Kn^ = 1, the curves 
for lift are similar, although slightly higher in magnitude, 
than those for drag. But as flow rarefaction increases, sen- 
sitivity to T w increases, and the distinction between the 
curves at each M^ begins to diminish. Indeed, these simu- 
lation results agree with free molecular flow theory 6 which 
predicts the following dependence, 


?— 2 


erf (J3S) + yj ttT^/T^PS 


Aerodynamic Sensitivity to Surface Temperature 

For a simple assessment of aerodynamic sensitivity to 
surface temperatures, steady-state rarefied flow past a flat 
thin plate was simulated with the DSMC technique using a 
computationally efficient code. 5 The two-dimensional sim- 
ulation geometry as used throughout the present work is de- 
picted in Fig. 1. The pitch angle was fixed at 0 = 45° and 
the intermolecular potential in the diatomic gas 3 was de- 
scribed by a = 8. The rotational collision number was set 
to unity while that for vibration was set to zero in order to 
maintain a fixed number of internal degrees of freedom of 
two and a constant ratio of specific heats, j = 1.4. Flow 
conditions were varied as follows, 

Kn^* {1,10,00} 


As demonstrated here, heating, drag, and especially lift 
are sensitive to surface temperature in hypersonic rarefied 
flows. Use of inaccurate estimates of T w in simulations of 
these flows could lead to unreliable aerodynamic predic- 
tions. 

Model for Dynamic Surface Temperatures 

For rarefied entry flows about thin panels, a new model 
is employed in the present work which couples the surface 
heat transfer into the flow solution to compute dynamic sur- 
face temperatures directly rather than prescribe some esti- 
mated constant temperature profile to the vehicle. A panel 
surface is composed of individual planar facets, each mod- 
eled with a thermal node on the front and back sides. Each 
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node i has an associated heat capacity c*, temperature 
radiating emissivity and dynamic net convective heat 
flux qi(t) per unit area; and is connected to neighboring 
nodes j with conductivity Kij. The corresponding energy 
balance for each node i is given by 

*(*)-«' (if - T^) ^ (T t - Tj) = (1) 

J 

where T D is the deep-space temperature. Coefficients c 
and k are written in a form which already accounts for ma- 
terial thicknesses and densities. 

For practical application to thin panels in the present 
work, we replace the general "z, j” network in Eq. (1) 
with independent dual-node surface facets. That is, ther- 
mal conductivity is considered only through the panel from 
front (side "a”) to back (side “6”) of a given facet, and is ne- 
glected laterally along the panel surface from facet to facet. 
Furthermore, at all but the smallest incidence angles j3 for 
panels in rarefied flow, convective heat flux on the back side 
may be neglected ( q b — ► 0). These assumptions lead to the 
following set of coupled equations which represent a dual- 
node model for the thermal characteristics of a thin panel 
surface facet in hypersonic rarefied flows, 

«.(0 - (Tt -K)~ k - T ‘> = C °^7T 

dT (2) 

-e h a (Tt -T* d ) + k (T a -T b ) = c b -^. 

This model may be applied regardless of the assumed ther- 
mal accommodation coefficient or surface reflection model 
(diffuse versus specular, etc.). Furthermore, a unique set of 
equations (2) may be written for each surface facet inde- 
pendently. Note that q represents the net convective energy 
flux, accounting for both incident and reflected molecular 
energies. The peak value Q for the net heat flux is simply 
measured directly in the DSMC code for each facet, and 
it accounts for any contributions of both translational and 
internal energy modes to gas-surface interaction. 

When considering Eqs. (1) or (2), it is important to 
distinguish time t during the aerobraking maneuver from 
the time-steps employed in a DSMC method. Note that 
the DSMC method typically simulates flow at a single in- 
stant in time only; say, at the perigee point, r, where peak 
heating occurs. Although the method will likely require 
several averaging-steps to yield a statistically meaningful 
solution, it nonetheless simulates the flow at a single in- 
stant during aerobraking. Equations (2) are coupled into 
the DSMC code by first assigning some initial tempera- 
ture profile over the vehicle surface, and measuring the net 
energy flux Q for each facet after one step. Then, by as- 
suming a time-dependent form for the heating pulse, q(t ), 


with peak value Q, Eqs. (2) may be integrated over t dur- 
ing aerobraking to solve for the dynamic surface temper- 
atures T a (t) and T a (t ), and thus for the temperatures at 
perigee, T a (r) and Ti(r). The front-side temperature T a (r) 
is then reassigned to the surface facet for use in the next 
DSMC step. After several steps, the coupling of this time- 
interation procedure with the DSMC method improves the 
estimates for Q and T a (r) simultaneously. 

Fortunately, the computational penalty associated with 
this procedure is negligible in comparison to computation 
of particle motion and interaction on each simulation step. 
Again, this procedure is not limited to thin panels, but may 
be applied readily to single-sided body surfaces. 4 

Application of Coupled Models 

Certain parameters in Eqs. (2) may dominate others to 
yield unique surface thermal characteristics. For example, 
if the heat capacities are extremely large, then the surfaces 
are effectively isothermal; such was the condition under 
which the temperature-sensitivity studies were performed 
above. If heat capacities are extremely small, then the panel 
becomes a conducting, double-sided radiator with instant 
thermal response to the dynamic heat flux. This condition 
still couples the surface thermal model to the flow-field sim- 
ulation, and will prove helpful in validating these models 
through comparison of free molecular predictions to theory. 

Perfect Conductors in Radiative Equilibrium 

If the heat capacity of the material is negligible (c — ► 0), 
then the panel exhibits instant thermal response to the heat 
flux. Consequently, surface temperatures for each facet 
may be assessed by replacing the transient flux q(t) by 
the steady flux measured in the DSMC code, Q . Further- 
more, when the thermal conductivity is large ( k —* oc), 
then surface temperatures on both sides of the panel must 
be equal (T a = T&). The surface facet therefore comes a 
double-sided node in radiative equilibrium with deep space, 
and Eqs. (2) reduce to a single analytic equation given by 

Q = (£a+e b )<T(Tt-T t D ). (3) 

Note that surfaces with zero thermal conductivity (perfect 
insulators) may be modeled by neglecting the radiative 
emission of the back side (setting e b = 0). 

To apply the model to a DSMC code, one must first as- 
sign some initial surface temperature to each surface facet. 
During a given simulation step, sample the net convective 
heat flux Q for each facet and update T a by solving Eq. (3) 
for each. This procedure is stable and introduces negligible 
additional computational burden even when computing T a 
for every surface facet after each DSMC step. More im- 
portantly, this procedure computes the temperature for each 
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facet independently and yields accurate surface temperature 
profiles which are consistent with the assumptions inherent 
to Eq. (3). This model was used in simulating aerobraking 
of the Magellan spacecraft 7 A similar, though considerably 
more cumbersome model has been suggested by Carey. 8 

To validate this coupled model, steady free molecular 
flow was simulated with the DSMC technique about flat 
panels at several flow conditions and angles of incidence. 
Simulation predictions of surface temperature were com- 
pared to those from free molecule theory 9 for perfect con- 
ductors in radiative equilibrium. Note that surface heating 
and temperature are uniform across the entire planar panel 
in free molecule flows. Simulation conditions are described 
by all combinations of the following parameters. 


dE a /dT a = -4 e a <rTl - k (8) 

dEi/dTb = -Ae b aTl - k. (9) 

The solution of Eq. (5) for one iteration leads to improved 
temperature estimates T' as follows, 

= ^ - g„ ( m/dTb ) 

° {dEa/dTa) ( dEb/dTb ) - k 2 

c T _ Eb + liST a 

*’ ( dEb/dTb ) 

T' a =Ta+STa 

Tl = Tb+6Tb. 


( 10 ) 

( 11 ) 

( 12 ) 

(13) 


Ea.b = 0.80 

0 = {5°, 15°, 30°, 45°, 60°, 90°} A =1.0 
p 00 = 3.72 x 10 -8 kg/m 3 7=1.4 

These cases were repeated using a thermal accommodation 
coefficient of A = 0.5 at Mach number M^ = 15. Re- 
sults are plotted in Fig. 5 and show excellent agreement for 
all cases. The following similarity parameter was used to 
normalize the surface temperature 9 



These results verify that the thermal model not only com- 
putes the correct temperature, but also that the temperature 
is coupled appropriately to yield the correct heating rate. 

Finite Conductors in Radiative Equilibrium 

Returning to Eqs. (2) with finite thermal conductivity 
through the panel, but maintaining zero heat capacities, 
yields a set of coupled quartic polynomials which may be 
solved for T a and T b for each double-sided surface facet 
after each DSMC step. Again, the instant thermal re- 
sponse of the material permits replacing q(t) by Q mea- 
sured in the DSMC code for each facet. Solution of the 
coupled polynomials is accomplished most easily through 
two-variable Newton iteration which involves solution of 
the vector equation, 

dE 

n 6T m “ E (5) 

where the error residuals E and their derivatives are defined 
from Eq. (2) as follows, 

Ea^Q-CaViTt-^J-KiTa-Tb) (6) 

Eb = -Sb<T {Tt -T * d )+k (Ta - Tb) (7) 


The iteration procedure involves using initial tempera- 
ture estimates T a and T b to evaluate E ai E by dE a /dT a > 
and dEb/dT bf and then computing improved estimates T f a 
and r/. Although this process could be repeated many 
times per DSMC time-step to yield accurate solutions, sev- 
eral such iterations are not necessary. Rather, performing 
only one Newton iteration step during each DSMC step will 
suffice; successive DSMC steps will rapidly drive the sur- 
face temperatures to their correct values. That is, as the 
DSMC solution progresses, the statistical results for Q im- 
prove as do the estimates for T a and 7i per facet. 

This Newton iteration process is demonstrated in Fig. 6 
which plots the evolution of T a and T b during the simulation 
for free molecule flow of Ni past a radiating, conducting, 
panel with conditions specified as follows, 

M oo =20 k = 20 W/m 2 K 

7 = 1.4 0 = 45° 

p x = 8 x 10 17 kg/m 3 A = 1.0 

T d = OK e a ,b = 0.80 

The initial guess for each temperature was the free-stream 
value, = 200 K. Despite jumping to nearly 117^ on the 
first iteration, the method quickly relaxed the temperatures 
down to their correct values after just ten DSMC steps. In- 
deed, this procedure is robust and stable as long as the initial 
guess for surface temperature is within roughly 2 orders of 
magnitude of the target value. 

Note in the figure that the difference between T a and T b 
is a result of finite thermal conductivity through the panel. 
This DSMC simulation was repeated using all combina- 
tions of the following parameters, 

M m = {5,10,20,30} 

/?= {15°, 45°, 75°} 
k = {2,20}W/m 2 K 
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Results of these DSMC simulations are plotted in Fig. 7 and 
compared to theory 9 using similarity parameters defined by 


H a - 





xA( noo kTj 3 ' 2 Z 


where 


Va 


ui = 5/2 + S 2 - 
S sin/? 


"exp (-S 2 sin 2 /?) 




4y/rZ 
1 + erf (5 sin /?) 


"exp (-S 2 sin 2 /3) 
2y/x 


(14) 

05) 


(16) 


(17) 


profiles obtained through numerical integration of Eqs. (2) 
for the same peak heating condition. The contribution from 
each mode of heat transfer is plotted. Note that the coupled 
DSMC temperature predictions compare well with the in- 
tegrated results. These temporal profiles demonstrate that 
the non-zero heat capacity of the panel material caused the 
peak temperatures to occur after peak heating. 

This simulation was repeated for = {5, 10, 15, 30} 
and the results were plotted in Fig. 9. DSMC surface tem- 
peratures compare well to analytic solutions from Eqs (2). 
Note from the analytic curves that, except at the highest and 
lowest energies, these temperatures are well below those 
corresponding to perfect insulators (k=0) or perfect con- 
ductors (/c=oc) in radiative equilibrium. These results not 
only validate the DSMC implementation of the dynamic, 
coupled, model from Eqs. (2), but also demonstrate the sig- 
nificant differences between the various models. 

Concluding Remarks 


Again, DSMC results for T a and 7* agree with theory over 
the entire flow range, validating the coupling of the DSMC 
method and the surface temperature model. 


Heat Capacity and Dynamic Thermal Response 

Non-zero material heat capacity leads to surface temper- 
atures which will respond dynamically to a transient con- 
vective heating pulse, q(t). Equations (2) were coupled into 
the DSMC technique to simulate free molecular flow of Nj 
past a thin panel under the following conditions. 


M co = 14 
7= 1.4 

= 5.00 x 10“ 7 kg/m 3 
T TO = 200 K 
j3 = 45° 

A = 1.0 


£a,b = 0.80 

k = 60 W/m 2 K 
c a>b = 2000.0 J/m 2 K 
r = 250 s 
T d = 150 K 


The form of the heating pulse was an empirical fit to an 
aerobraking profile 7 and is given by 

q(t) = “Q ^1 - cos exp ^—2 - 2cos t (18) 

where r represents a time constant such that peak heating 
occurs at i = r and the duration of the entire heating pulse 
is 2 r. The initial surface temperatures, T a ( 0) and 7i(0), are 
assumed to equal the deep-space temperature T D . 

The simulation predicted a peak heating value of Q = 
11529.4 W/m 2 at perigee and a front-side perigee temper- 
ature of T a = 370.8 K and a peak temperature of T a = 
460.5 K during the heating pulse. These points are plotted 
in Fig. 8 along with the temporal heating and temperature 


These results demonstrate that thermal conductivity and 
heat capacity play significant roles in establishing the sur- 
face temperatures of thin panels in rarefied hypersonic flow 
during dynamic aeropass maneuvers. DSMC simulations 
of flows past isothermal plates proved that under many con- 
ditions, particularly in highly-rarefied flows, surface tem- 
perature has great impact upon vehicle aerodynamics and 
heat transfer. It follows that accurate simulation of aero- 
braking scenarios requires accurate assessment of these sur- 
face temperatures. For dynamic heating pulses, this is af- 
forded by time-integration of the surface heat transfer equa- 
tions (2) and coupling of the surface temperatures back into 
the flow simulation. As such, surface temperatures are com- 
puted directly in the DSMC technique rather than relying 
upon isothermal estimates. These models may be employed 
in continuum methods for computational fluid mechanics 
as well, although the effects of surface temperature may be 
less significant in the continuum flow regime. 

The differential equations could be simpified for mate- 
rials with zero heat capacity, representing a double-sided 
panel with finite thermal conductivity from front to back. 
This coupled model was validated through comparison of 
simulated surface temperatures over a large range of flow 
conditions to free molecule theory under the same condi- 
tions. Likewise, the models may be simplified further for 
perfect conductors in radiative equilibrium. Again, com- 
parisons with free molecule theory validated this model 
over a large range of flow conditions. Naturally, these 
DSMC models may be applied regardless of Kn^ ; free 
molecule flow was studied here for validation purposes due 
to the ready availability of the corresponding analytic solu- 
tions. 
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Fig. 1 Schematic represen tatio a of flow simulation geometry (not 
to scale). 
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Fig. 2 Effects of M w , and Kn^ upon net convective heat- 
ing?. 
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Fig. 3 Effects of T w , , and Kn^ upon drag coefficient. 





Fig. 4 Effects of T w , M m p and upon lift coefficient 
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Fig- 5 Results from DSMC and theory for radiating, perfectly- 
conducting panels in free molecule flow. 



Fig. 6 Relaxation history for surface temperature computed in the 
DSMC code with the single-Newton -iteration scheme for radiat- 
ing, conducting panels. 



Fig. 7 Results from DSMC and theory for radiating, finitely- 
conducting panels in free molecule flow, 




Fig. 8 Transient results from Eq. (2) for: (a) heating due to con- 
vection, conduction, and radiation; and (b) transient temperatures 
on each side of the panel. Temperatures are compared to DSMC 
results employing this model. 



Fig. 9 Perigee and peak surface temperatures for front and back 
sides of a panel over a large range of peak convective heating. So- 
lutions of Eq. (2) are compared to DSMC results using this model 
and to perfect conductors and insulators in radiative equilibrium. 





Appendix C 





AIAA 93-2806 

Flow Resolution and Domain of Influence 
in Rarefied Hypersonic Blunt-body Flows 

Brian L. Haas 
Eloret Institute 
Palo Alto, CA 94303 


AIAA 28th Thermophysics Conference 

July 6-9, 1 993 / Orlando, FL 


For permission to copy or republish, contact the American Institute of Aeronautics and Astronautics 
370 L'Enfant Promenade, S.W., Washington, D.C. 20024 



AJAA 93-2806 1 


FLOW RESOLUTION AND DOMAIN OF INFLUENCE 
IN RAREFIED HYPERSONIC BLUNT-BODY FLOWS 

Brian L. Haas* 

Eloret Institute , Palo Alto, California 94303 

This parametric study assesses the influence of upstream domain size and grid resolution upon flow 
properties and body aerodynamics computed for rarefied flows over cold blunt bodies with a direct sim- 
ulation Monte Carlo (DSMC) particle method. Empirical correlations are suggested for aerodynamic co- 
efficients for two-dimensional flows past a perpendicular flat plate. Free-stream parameters which were 
varied in the study include the Mach number, Knudsen number, surface temperature, and intermolecular 
potential. Insufficient grid resolution leads to overprediction of aerodynamic heating and forces in the 
DSMC method. Solution accuracy correlates well with the Reynolds number defined at the wall temper- 
ature and the stagnation mean free path relative to the cell dimension. Molecular diffusion firm the body 
surface results in greater far- field domain influence as flow rare faction increases. As a result, insufficient 
upstream domain size in the DSMC method leads to overprediction of heating and drag. Errors in aerody- 
namic coefficients correlate well with x, the distance ahead of the body where flow temperature reaches 
half of its peak value. Engineering criteria for DSMC cell resolution and domain size are suggested. 
Simulation of a hard-sphere gas is more sensitive to grid resolution while simulation of a Maxwell gas is 
more sensitive to upstream domain size. 


NOMENCLATURE 

C aerodynamic force coefficient 

H plate height 

Kn Knudsen number, X/H 

Kn Whitfield’s correlation parameter, Eq. (6) 

L upstream domain size 

T domain size correlation parameter, Eq. (18) 

M Mach number 

Nu w Nusselt number at the wall, Eq. (7) 
n gas number density 

Pr Prandtl number, Eq. (8) 

Q free-stream kinetic heat flux, u 3 ^ /2 

q aerodynamic net heat flux 

Re Reynolds number, Eqs. (2) and (3) 

St Stanton number, Eq. (9) 

T temperature 

u flow velocity 

x location relative to stagnation point 

x distance from plate to T half-rise location 

a intermolecular potential exponent 

7 ratio of specific heats = 1 .4 

6 length of a cubic grid cell 

A gas mean free path 

p gas mass density 

Superscripts: 

* presumed-correct value 

DOM pertains to domain effects 

res pertains to resolution effects 


* Research Scientist Member, AIAA. 

Mailing Address: NASA Ames Research Center 
M/5 230-2, Moffett Held, California 94035-1000 

This piper is declared a work of the US. Government and is not 
subject to copyright protection in the United States. 


Subscripts: 

a 

adiabatic wall condition 

D 

drag 

ERR 

error coefficient 

L 

lift over top half of body 

s 

stagnation value 

w 

wall value 

oo 

freestream value 


INTRODUCTION 

Row field characteristics about blunt bodies during at- 
mospheric entry lead to considerable challenges for com- 
putational simulation. Highly rarefied flows, where the gas 
density n is low and the mean free path A between molecu- 
lar collisions is high, are better suited to particle simula- 
tion methods, such as the direct simulation Monte Carlo 
(DSMC) technique pioneered by Bird, 1 rather than contin- 
uum techniques based upon the Navier-Stokes equations. 
DSMC methods employ many model particles whose mo- 
tion and interaction simulate gas dynamics directly. The 
simulated flowfield is divided into a network of small cells 
to facilitate collision modeling and statistical sampling. 

The computational burden of DSMC methods, however, 
grows proportionally with local gas density and the size 
of the computational domain. In typical entry flows, the 
body surface temperature is rather low compared with the 
stagnation temperature, leading to a steep density gradient 
near the body surface. Accurate simulation of this flow 
requires sufficient grid resolution near the body. Alterna- 
tively, the extent of freestream rarefaction results in a lead- 
ing shock layer which is fully merged with the boundary 
layer of the body, yet extends far upstream. Accurate sim- 
ulation therefore requires a large upstream computational 
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flow domain. Together, these requirements rapidly drive 
the computational cost of the particle simulation method 
upward. Hermina 2 and Woronowicz 3 studied grid resolu- 
tion effects for flows over slender bodies and flat plates par- 
allel to the free-stream direction. In contrast, the present 
study addresses flow over blunt bodies. Also, errors result- 
ing from limited domain size became evident in recent com- 
putations of aerobraking of the Magellan spacecraft 4 . The 
DSMC researcher must therefore understand what simula- 
tion costs are necessary to obtain a solution of sufficient 
accuracy. 

The objective of the present parameteric study is to as- 
sess quantitatively the sensitivity of aerodynamic loads and 
gas properties to grid resolution and simulation domain size 
in blunt-body, hypersonic, rarefied flows. In addition, cor- 
relations are suggested for aerodynamic coefficients and 
flow properties over a wide range of free-stream conditions. 

METHODS OF ANALYSIS 

The present study employed an efficient particle simula- 
tion technique 5,6 to investigate the net heat flux q and aero- 
dynamic forces on a two-dimensional flat plate of height H 
facing perpendicular to the flow direction. Gas tempera- 
ture T and number density n were also computed along the 
stagnation streamline. The following range of flow condi- 
tions were used: 


M oo = {5,10,20} 

Kn TO = {0.1, 0.3, 1,3, 10} 

T tt /T 00 = { 1,5,9} 
a = {4, 8 , oo}. 

For each flow case, the gas was modeled with two fixed in- 
tonal degrees of freedom leading to a constant ratio of spe- 
cific heats of 7 = 1.4. Interaction of the gas with the plate 
surface at wall temperature T w was modeled as fully diffuse 
and thermally accommodated. All simulations employed 
between 4 and 64 particles per cell as the input freestream 
number density, . The code used cubic cartesian cells of 
uniform dimension throughout the flowfield, but this does 
not alter the conclusions of the present work when applied 
to DSMC codes with variable cell dimensions. 

For each set of flow conditions, the simulation resolution 
and domain were varied considerably and their effects upon 
solution accuracy were determined. Grid resolution is de- 
fined by the plate height H while the flow domain size is de- 
fined by the distance L of the computational field upstream 
of the plate, where both H and L are measured in units of 
grid cells of linear dimension 6. These length scales are 
depicted schematically in Fig. 1 along with x which repre- 
sents the location, far upstream of the plate, where the local 
translational temperature first reaches half of its peak value. 

Estimate of Local Mean Free Path 

In the DSMC technique, two particles in the flow may 
collide only if they both reside in the same grid cell. Fur- 


thermore, any two particles within a given cell can collide 
regardless of their respective positions in the cell. The colli- 
sion probability for a particular pair of neighboring particles 
is assessed as a function of molecular constants and the rela- 
tive translational energy of the pair . 5 If the cell dimension <5 
near the cold body surface is too large, then energetic par- 
ticles at the far edge of the cell are able to readily transmit 
momentum and energy to particles immediately adjacent to 
the surface. The latter particles may, in turn, transmit that 
energy and momentum to the surface. This leads to over- 
prediction of both the surface heat flux and the aerodynamic 
force on the body than would occur in the real gas. This 
error is minimized by reducing the cell dimension relative 
to the local mean free path of molecules near the surface. 
From equilibrium kinetic theory, where translational veloc- 
ities have a Maxwellian distribution, the local mean free 
path A is related to the local number density n and temper- 
ature T as follows. 


-fir) (0 


2/a 


( 1 ) 


The exponent of the assumed inverse-power intermolecu- 
lar potential, a, may vary between the limits of the Maxwell 
molecule (a = 4) and the hard sphere (a = oo). The expres- 
sion in Eq. (1) is not exact for flows far from translational 
equilibrium, but will suffice to characterize the degree of 
rarefaction in the flow. Employing the stagnation tempera- 
ture and density near the plate surface from the DSMC sim- 
ulation leads to an approximation of the stagnation mean 
free path, A,. Theoretical estimates of A, were derived 
by Whitfield , 7 but were limited to particular flow condi- 
tions and were not as convenient as the simple expression 
in Eq. (1) for use in the present study. 

As a consequence of Eq. (I), regions ofhigh density lead 
to short mean free paths, requiring finer cell resolution to 
capture flow gradients accurately. Such behavior charac- 
terizes hypersonic rarefied flows near cold blunt-body sur- 
faces. 


FLOW CHARACTERISTICS AND CORRELATIONS 

Having computed flow properties and aerodynamics 
over a wide range of input conditions and simulation param- 
eters, it proved instructive to summarize the most accurate 
results; that is, those corresponding to grids of sufficient 
resolution and domains of sufficient size. 


Flow Along the Stagnation Streamline 

Density profiles along the stagnation streamline up- 
stream of the plate are plotted in Fig. 2a to demonstrate 
the effects of free-stream rarefaction. Note that the den- 
sity at the plate surface (at x/H = 0) is significantly 
higher than the free-stream value. This behavior is slightly 
more pronounced, and the gradients are steeper, for flows at 
lower Kn^ . This, combined with the fact that lower Kn^ 
means that the freestream A^ is short, dictates that the cell 
resolution must be very fine compared with more rarefied 
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flows. Notice that there is no real shock structure ahead of 
the body except at the lowest Kn^ . 

Temperature profiles along the stagnation streamline for 
these same flows are plotted in Fig. 2b. Temperature was 
defined by the variance of particle velocities. These curves 
clearly demonstrate that the body influences the flow far 
upstream. This domain of influence increases with Kn^ as 
a result of upstream diffusion of panicles which reflected 
from the body surface. By contrast, Kn^ has little effect 
upon the resultant stagnation temperature and density at the 
plate surface. 

For a given free-stream Knudsen number, the effects of 
Mach number on the streamline profiles is demonstrated in 
Fig. 3. Greater Mach number leads to significantly greater 
stagnation density and temperature, as well as greater peak 
temperature. However, unlike Kn^ , M^ has very little 
effect on the extent of the domain of influence upstream of 
the plate. 

As is evident in the profiles in Fig. 4, however, wall tem- 
perature T w affects all flow properties including the peak 
temperatures, stagnation values, and the extent of the do- 
main of influence of plate in these rarefied flows. This 
results from particles reflecting from hotter surfaces with 
greater energies, leading to greater relative translational ve- 
locities and flow temperatures. Energetic reflections from 
the hotter surface reduces the net buildup of particle density 
near the plate. 

Finally, the effects of employing different exponents a 
of the intermolecular potential are demonstrated in Fig. 5. 
Just as with Kn^, different potentials have little effect 
upon the stagnation temperature or density, but lower val- 
ues of a greatly increase the extent of the domain of in- 
fluence upstream of the body. Plotting the local mean free 
path as computed with Eq. (1) shows that A is much greater 
for a = 4 than for a = oo. The flow therefore becomes 
effectively more rarefied, and the profiles are indicative of 
those corresponding to large Kn^ in Fig. 2. Also, the gra- 
dients for the hard-sphere are steeper and lead to lower A, at 
the stagnation point on the body compared to the Maxwell- 
molecule. As a consequence, accurate DSMC solutions for 
hard-sphere gases require finer cell resolution while solu- 
tions for Maxwell-molecule gases require larger upstream 
computational fields. 

Many of these effects are summarized by plotting A, and 
x for all flow conditions. The stagnation mean free path 
correlated fairly well with the product Re^Kn^ 4 as plotted 
in Fig. 6. The Reynolds number based on the wall tem- 
perature is defined from the temperature dependence of the 
coefficient of viscosity and is given by 


Re^ = Re 


fZk 

V T. 


1 / 2 — 2 / a 


( 2 ) 


where the free-stream Reynolds number is determined from 
Bird’s definition of mean free path 8 as follows, 
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This correlation takes into account the effects of all flow pa- 
rameters (M to , Kn^ , a, and T w ) and demonstrates that A, 
is strongly dependent upon and weakly upon Kn^ . 
Employing a least-squares log-log fit to the data in Fig. 6 
leads to the following simple expression for the stagnation 
mean free path for these flows, 

\ - ^2° ('4') 

* 2Re U) Kn^/ 4 1 ’ 

As plotted in Fig. 7, the upstream extent x of the domain 
of influence of the body is correlated against the product, 

f_L_VM ‘ /4 (5) 

\l-2/aJ\T „) Re., 

Again, this gives a fair correlation and encompasses all pa- 
rameters. Here, x exhibits strong dependence upon Kn^ 
than upon M^ . 

Ail correlations suggested in the present study were de- 
veloped through parameterization. That is, simulation re- 
sults for cases with identical flow conditions except for a 
given parameter were compared with each other to assess 
dependence upon the given parameter. Simulation results 
were plotted using quantities multiplied by powers of the 
independent parameters. 


Correlations for Aerodynamic Coefficients 

No existing theory describes completely the heating and 
aerodynamic loading of bodies over the entire range of 
Reynolds number. The objective here is only to estimate 
the dependencies of the aerodynamic coefficients upon flow 
conditions empirically, rather than develop a detailed the- 
ory for them. In addition to numerous experimental investi- 
gations, previous computational attempts to identify corre- 
lation parameters have met with some success, including 
those of Woronowicz 9 for flat plates parallel to the free- 
stream flow, and Gilmore and Harvey 10 for flat-ended cylin- 
ders. The latter work concludes that Whitfield’s parame- 
ter Kn, defined by 

ra '^~(\/^l lM “ +i ) '■ <6> 

is reasonable for correlating heat transfer. The results of 
the present study employed this parameter to correlate the 
wall Nusselt number, NUu,, divided by ^/Kn^, to obtain 
the result plotted in Fig. 8. Here, Niv, is defined from the 
following expressions 


NUu, = Re^ Pr St 
47 


Pr = 


St = 


97-5 
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q/Q 
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The last expression closely approximates the adiabatic wall 
temperature for free molecule flows 11 and is assumed to be 
a fair estimate for the rarefied flows studied here. 

Aerodynamic drag and lift on the plate did not vary sig- 
nificantly among the flow cases, but did exhibit some weak 
correlation to the ratio Re^ /Kn^ as plotted in Fig. 9. In 
the present context, the lift coefficient C L is defined by the 
body force, perpendicular to the flow direction, over the top 
half of the plate only. Defining C L in this way yields a pa- 
rameter which will prove below to be sensitive to grid reso- 
lution and domain size. Naturally, integrating the force over 
the top and bottom halves of the plate would yield zero lift 
force due to symmetry. 

GRID RESOLUTION 

For all flow conditions within the Knudsen number 
range, 0. 1 < Kn^ < 1 .0, simulations were performed 
to assess the effects of grid resolution upon solution ac- 
curacy. Since all flow parameters and aerodynamic co- 
efficients are normalized by the plate height, grid resolu- 
tion was defined by the ratio H/8. At each flow con- 
dition, simulations were performed using plates measur- 
ing H/8 = {2, 4, 6, 10, 20, 40, 60, 100, 200, 400} . 

Resolution Effects On Flow Profiles 

The density and temperature profiles along the stagna- 
tion streamline ahead of the body are plotted in Fig. 10. 
For the coarse grids, where H/8< 100, the profiles were 
smeared considerably. For H/8> 100, the profiles coa- 
lesced to a single form and were assumed therefore to rep- 
resent accurate solutions. This same general behavior was 
observed for all flow conditions, although the resolution 
required to capture the appropriate profiles did not need 
to be as fine at larger Kn TO compared with lower Kn^ . 
For example, at Kn^ =1.0, accurate profiles were obtained 
at H/8> 60. 

Resolution Effects On Aerodynamic Coefficients 

For each simulation, the coefficients for heating and 
aerodynamic forces were computed. To assess the errors 
corresponding to a given flow resolution, the aerodynamic 
coefficients computed with the finest resolution at each flow 
condition were assumed to represent the “correct” values 
for that flow. It was the stagnation density and temperature 
from that solution which was also used to determine A, as 
plotted in Fig. 6. 

As observed qualitatively for all flow conditions, coarse 
grid resolution led to overprediction of heating, drag, and 
lift The errors dropped considerably with finer resolution. 
Errors in the aerodynamic coefficients were correlated and 
plotted against the product A, Re^ /8 in Fig. 1 1 . Best agree- 
ment was obtained, through parameterization, by employ- 
ing error coefficients defined as follows. 


(Cd) „, . 5 «/. (<* _ , j (12) 

(C.C-IC 1 ' 5 (§*-■)• (13) 

Superscript denotes the value for the coefficient which 
was assumed to be most accurate from the finest resolution 
simulation. These plots provide an engineering assessment 
of the accuracy associated with a given grid resolution, in 
terms of stagnation mean free path, for a given set of flow 
conditions. Scatter in the plots resulted from imperfect cor- 
relation parameters and the statistical limitations inherent 
to the DSMC technique. Scatter is particularly noticeable 
at the highest resolutions where the errors are small relative 
to the accuracy of the simulation method. 

A suitable criterion for engineering applications, which 
estimates that the error coefficients will drop by two orders 
of magnitude, specifies a cell resolution near the body sur- 
face dictated by the stagnation mean free path as follows, 

Re 

6 < \, ~ 2 &- Resolution Criterion (14) 

Employing the linear fit from Eq. (4) leads to the criterion 
in teims of flow input parameters only, 

6 < 4 q^ 3/4 • Resolution Criterion (15) 

*00 

DOMAIN SIZE 

For rarefied flows in the higher Knudsen number 
range (l.O^Kn^lO.O), simulations were performed to 
assess the influence of the upstream flow field domain upon 
solution accuracy. In view of the results above, the cell 
resolution used in each case here was sufficiently fine such 
that the stagnation mean free path exceeded one cell -length 
near the body. The objective here was to employ different 
upstream computational field sizes, defined by length L in 
Fig. 1, to determine the minimum acceptable domain size 
required for a given penalty in solution accuracy. At each 
flow condition, simulations were performed using domains 
of dimensions L/H = {0.1, 0.2, 0.5, 1,2,4, 8, 16,32}. 

Domain Size Effects On Flow Profiles 

Flow density and temperature along the stagnation 
streamline are plotted for one flow condition in Fig. 12 em- 
ploying different domain sizes L. With sufficiently large L, 
all the curves coalesced to a single curve to represent the 
assumed correct profile as represented at L/H = 32. With 
insufficient upstream domain (i.e. small L), density and up- 
stream temperature are underpredicted while the peak tem- 
perature is overpredicted. However, once the simulation 
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domain size exceeded roughly twice the distance of the tem- 
perature half-rise point, x, the temperature and density pro- 
files followed the expected forms predicted when ample 
upstream domain was employed. Although not presented 
here, behavior for all other flow conditions are qualitatively 
similar to these results. 


Domain Size Effects On Aerodynamic Coefficients 

Quantitative assessment of the errors in aerodynamic 
loads on the plate for simulations employing different up- 
stream domain sizes are presented in Fig. 13. The errors for 
each case are computed relative to the aerodyanmic loads 
which were obtained using the largest domain L and pre- 
sumably represent the correct values. Each plot employs 
a length scale normalized by x. As expected, the errors in 
predicted aerodynamics decreased as larger upstream do- 
main sizes were used. Aerodynamic heating q and drag C 0 
were overpredicted in all cases. The lift forces were so 
small for these rarefied flows that the errors for each domain 
size were negligible relative to the statistical fluctuation of 
the simulation, and were therefore not plotted. Through pa- 
rameterization, the best correlations over the range of flow 
conditions employed error coefficients defined by. 



and a domain-size correlation -parameter given by 



Note that the data levels off as the domain size parame- 
ter L drops below roughly 0.3. In the limit of decreasing up- 
stream domain, the incoming flow is almost completely un- 
affected by particles reflecting from the plate surface such 
that the aerodynamic loads approach the values correspond- 
ing to the free molecular flow limit. For heating and drag, 
the limit is less erroneous for flows at greater frees tream 
Knudsen numbers Kn^ since free molecular flow repre- 
sents the limit of greatest possible rarefaction. The inter- 
esting paradox is that, while the domain of upstream in- 
fluence increases with Kn^ , the penalty in terms of solu- 
tion accuracy is less dramatic compared to simulations at 
lower Kn^ . 

For engineering purposes, an appropriate criterion for 
defining a suitable computational flow domain for rarefied 
flow simulation is given when the error coefficients drop by 
two orders of magnitude, leading to the criterion, 


L > 2. 


CONCLUDING REMARKS 

When employing direct particle simulation methods, ac- 
curate simulation of highly rarefied flows about cold blunt 
bodies requires sufficient computational grid resolution and 
upstream flow domain size. The objective of this paramet- 
ric study was to assess quantitatively what penalty is suf- 
fered in solution accuracy when minimizing computational 
expense by using short domains and coarse grids. Qual- 
itative results were similar for all free-stream conditions. 
Criteria for grid resolution and domain size were expressed 
in terms of flow properties A, and x. These criteria apply 
specifically to the two-dimensional flow about a flat plate 
perpendicular to the freestream direction when 7 = 1.4. 
However, these results likely apply to any cold- wall, blunt- 
body, hypersonic, rarefied flow problem. Empirical corre- 
lations were suggested for stagnation streamline properties 
and aerodynamic loading over the range of flow conditions. 
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Fig. 1 Definition of simulation length-scales and schematic 
representation of temperature and density profiles upstream 
of a cold plate of dimension H in rarefied flow. 
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Fig. 2 Effects of flow rarefaction, Kn^ , upon profiles along 
the stagnation streamline; (a) number density, n; (b) tem- 
perature, T. 



Fig. 3 Effects of Mach number, M^ , upon profiles along 
the stagnation streamline. 



Fig. 4 Effects of wall temperature, T w , upon profiles along 
the stagnation streamline. 
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Location Upstream of Plate, x / H 


Fig, 5 Effects of intermolecular potential, a, upon profiles 
along the stagnation streamline. 



Fig. 6 Empirical correlation for the stagnation mean free 
path. 



Fig. 7 Empirical correlation for the extent of the domain of 
influence upstream of the body. 



Kn 

Fig. 8 Empirical correlation for plate heating. 



Fte w / Kru 

Fig. 9 Empirical correlation for plate drag and lift forces. 



Location Upstream of Plate, x / H 


Fig. 10 Effect of grid resolution upon profiles along the 
stagnation streamline upstream of the plate. 
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Fig. 1 1 Correlated effect of DSMC grid resolution upon aerodynamic coefficients 
for all flow conditions; (a) heating; (b) drag; (c) lift 
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Fig. 13 Correlated effect of DSMC domain size upon aerodynamic coefficients for 
all flow conditions; (a) heating; (b) drag. 
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SIMULATED RAREFIED AERODYNAMICS 
OF THE MAGELLAN SPACECRAFT DURING AEROBRAKING 

Brian L. Haas* 

Eloret Institute , Palo Alto , California 94303 
Dunvin A. Schmittf 

Martin Marietta Astronautics Company , Denver, Colorado 

Aerodynamic loads upon the Magellan spacecraft during aerobraking through the atmosphere of Venus 
are computed at off- design attitudes with a direct simulation Monte Carlo (DSMC) particle method . This 
method is not restricted to the assumption of collisionless flow normally employed to assess spacecraft 
aerodynamics. Simulated rarefied flows at nominal altitudes near 140 km and an entry speed of 8.6 km/s 
were compared to simulated and analytic free molecular results. Aerodynamic moments , forces, and 
heating for rarefied entry at all attitudes were 7- 1 0% below free molecular results. All moments acted to 
restore the vehicle to its nominal zero-pitch, zero-yaw attitude. Suggested canting of the solar panels is an 
innovative configuration to assess gas-surface interaction during aerobraking. The resulting roll torques 
about the central body-axis as predicted in rarefied flow simulations were nearly twice that predicted for 
free molecular How, although differences became less distinct for thermal accommodation coefficients 
well below unity. In general, roll torques increased dramatically with reduced accommodation coeffi- 
cients employed in the simulation. In the DSMC code , periodic free-molecule boundary conditions and 
a coarse computational grid and body resolution served to minimize the simulation size and cost while 
retaining solution validity. 


NOMENCLATURE 

A j vehicle frontal area (236.8 cells 2 ) 

a thermal accommodation coefficient 

D diameter of the HGA (370 cm, 11.77 cells) 

F aerodynamic force coefficient 

Kn Knudsen number, X/D 

M mach number or moment coefficient 

m molar mass of atmospheric gas 

7i number density 

Q incident free-stream energy flux, p^ /2 

q net heat flux on solar panels 

T temperature 

u flow velocity 

x-axis roll axis, along central body axis 

y-axis yaw axis, through solar panel axes 

^-axis pitch axis 

a intermolecular potential exponent 

i surface radiative emissivity 

A gas mean free path 

p gas mass density (kg/m 3 ) 

oo free-stream value 

Geometric Components: 

ALTA altimeter antenna 

BUS equipment bus 

FEM forward equipment module 
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Mailing Address: NASA Ames Research Center 
M/S 230-2, Moffett Field, California 94035-1000 
f Retired. Currently Vice President, Tactical Technical Solutions 
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subject to copyright protection in the United Stales. 


HGA high gain antenna 

LGA low gain antenna 

MGA medium gain antenna 

PAN solar panels 

REM rocket engine module 

SRM housing from solid rocket motor 


INTRODUCTION 

The Magellan spacecraft has been mapping the surface 
of Venus from a highly elliptic orbit (eccentricity, e = 0.39) 
since September of 1990. Mission planners at NASA’s Jet 
Propulsion Laboratory (JPL) would like to circularize the 
orbit to improve mapping but cannot perform the manuever 
through thruster activity alone due to limited remaining pro- 
pellant. The orbital maneuver will be achieved through a 
series of gentle passes through the Venus atmosphere (ap- 
proximately 840 passes with an average velocity decrement 
of 1.5 m/s each) 1 which were initiated on May 27, 1993. 
Besides circularizing the orbit, these maneuvers will pro- 
vide considerable data related to atmospheric entry of satel- 
lite spacecraft. 

The spacecraft configuration is depicted in Fig. 1 in its 
nominal entry orientation (three-axis stabilized) with flow 
directed along the central body axis (x-axis) such that the 
high gain antenna is aft. The solar panels in this config- 
uration are normal to the flow, but may be canted at any 
angle (about y-axis) to form an effective “windmill” dur- 
ing entry. The restoring roll torques on the spacecraft may 
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be measured along with surface heating and orbit altitudes 
to deduce the flow-field density and the surface accommo- 
dation coefficients in the normal and tangenual directions. 
Axes for pitch and yaw are indicated in the front view. 

One critical mission constraint is the aerodynamic heat- 
ing on the delicate solar panels during each aeropass. 
Direct particle simulations of entry at several altitudes 
from 125 km to 140 km in the nominal orientation 2 verified 
that this heating was within the specified tolerance at alti- 
tudes exceeding 135 km. At the nominal altitude ofl40 km, 
the heating and drag results corresponding to the rarefied 
flow were very close to those corresponding to free molec- 
ular (collisionless) flow at that altitude. 

Mission planners are also concerned with spacecraft sta- 
bility during aerobraking, particularly if the entry orienta- 
tion has a high angle of pitch or yaw. Simplified aerody- 
namic analyses rely upon the assumption of free molecular 
flow and may be computed reliably with the FREEMAC 
code. 3 However, this assumption may not be valid at nom- 
inal flight conditions when the free-stream Knudsen num- 
ber is below ten (based upon the HGA diameter). The ob- 
jective of the present study is to calculate torques, forces, 
and heating of the Magellan spacecraft for several entry at- 
titudes employing a DSMC particle method which is not 
restricted to the assumption of free molecule flow. The re- 
sults of the rarefied simulation may then be compared to 
FREEMAC results to assess the significance of molecular 
collisions in the flow. A similar study by Rault 4 employed 
a DSMC method and a very fine body resolution. Unfor- 
tunately, some of those results are questionable since they 
differed significantly from FREEMAC results even when 
simulating free molecule conditions. 

Roll torques associated with the windmill entry configu- 
rations were also computed in the present work at nominal 
conditions and compared to free molecule results. Aero- 
dynamic sensitivity to the assumed surface thermal accom- 
modation coefficient was assessed for the windmill config- 
uration. The nominal altitudes were 136—140 km and the 
periapsis velocity was u ^ = 8.6 km/s. 

PARTICLE SIMULATION METHOD 

Direct simulation Monte Carlo (DSMC) particle meth- 
ods model the motion and interaction of thousands or mil- 
lions of computational particles to simulate gas dynamics. 5 
Specifics of the DSMC code employed in these studies are 
provided in Ref. 2. Given a particular position, velocity, 
and internal energy status, each particle in the flow-field 
travels unobstructed along the linear trajectory of its veloc- 
ity vector over the duration of a single ti me step . Neigh- 
boring particles are then identified throughout the flow-field 
and paired-off as potential collision candidates. The flow- 
field is divided into a network of uniform cubic cells to fa- 
cilitate identification of neighboring particles and to define 
the finest resolution for sampling macroscopic flow prop- 
erties. Employing probabilities as functions of individual 
collision parameters such as collision cross-section and rel- 


ative translational speed 6 , the subset of all candidate pairs 
which collide during the timestep are identified. In simulat- 
ing free molecule flow, the probability of collision is arti- 
ficially set to zero, representing an infinite molecular mean 
free path. 

The entire simulated flow-field is initialized with free 
stream conditions. The particle simulation then runs 
through a transient phase as the solution develops and flow- 
field structures form. Upon reaching steady-state, the sim- 
ulation collects statistical samples for measuring properties 
of the flow-field and body surface fluxes. 

The DSMC code employed in the present study was 
developed by McDonald 7 for efficient implementation on 
vector supercomputers. This code simulates non-reactive, 
three-dimensional, flow of general gas mixtures about arbi- 
trary geometries. Molecular collisions pertain to the Vari- 
able Hard Sphere model of Bird 8,6 with an inverse inter- 
molecular potential exponent of a = 5 to simulate the gas 
dominated by CO 2 . The internal energy modes are mod- 
eled with three fully-excited degrees of freedom to account 
for molecular rotation and vibration. Internal energy exci- 
tation is performed with the mechanics of Borgnakke and 
Larsen 9 employing a fixed probability of relaxation of 1/5. 

Body Geometry and Grid 

The geometry of the Magellan spacecraft is shown in 
Fig. 1 and compared to the simulated geometry. Due to 
the large mean free path in the flow about the spacecraft, 
small features on the vehicle such as the altimeter antenna 
(ALTA), the medium and low gain antennas (MG A, LGA), 
and the rocket engine modules (REM) have negligible im- 
pact on the flow-field as a whole and may be excluded to 
simplify the simulation geometry, leading to two planes of 
symmetry on the vehicle. These small features do impact 
the vehicle aerodynamics, however, and their contribution 
will be considered separately. 

To represent a surface in the cubic cartesian grid net- 
work of the simulation, it is necessary to approximate the 
surface as a composite of planar facets. Each facet has a 
normal defined from the intersections of the surface with 
the edges of the cell. This faceted description of the body 
is appropriate given that body radii are large in comparison 
to the cell size, and that the intersection of different sur- 
faces occurs at cell boundaries. Since the solar panels and 
HGA are dominant components of the structure, the sim- 
ulation cell network was scaled such that the square sim- 
ulated panels (measuring 8 cells x 8 cells) have the same 
frontal area as the actual solar panels and that the HGA di- 
ameter (11.77 cells) matches the actual HGA via scaling 
factor, 1 cell = 31.44 cm. Such a coarse grid and body res- 
olution was used in this study to reduce the computational 
expense associated with running several simulations. 

Grid resolution greatly impacts the accuracy of DSMC 
solutions for vehicle aerodynamics and heating in rarefied 
flows. An established criterion for sufficient grid resolu- 
tion is that the local mean free path must exceed the cell 
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dimension. For cold-wall blunt-body rarefied flows, flow- 
field density will likely be quite large near the body surface, 
leading to a small stagnation mean free path. Equilibrium 
kinetic theory provides the following estimate of the local 
mean free path 10 near the body surface, 


A oo n \ T oo J 


( 1 ) 


Simulated at nominal flight conditions, the flow density n, 
local mean free path A, and temperature T along the stag- 
nation streamline at the center of the solar panel are plotted 
in Figs. 2 and 3. Flow properties at the panel surface led to 
stagnation mean free paths in the range 1.8 < A < 4.2 cells, 
exceeding the minimum accuracy criterion at all altitudes. 

Body surfaces are modeled 11 ■ 12 as if in radiative equilib- 
rium with deep space at a temperature of 4 K and with an 
emissivity of c = 0.82. Advantages of this model include: 
1) it is simple to implement in the simulation; 2) it does not 
require a prescribed estimate of surface temperature; and 
3) it allows each surface facet to reach its own temperature 
independent of neighboring facets. It is assumed that radia- 
tion from the flow-field or from other body surfaces would 
contribute negligibly to the net heating of a given surface 
facet. Most results will employ a thermal accommodation 
coefficient of unity, typical of rough cool surfaces facing 
into the flow. Molecules reflect diffusively from the sur- 
face with translational and internal energies sampled at the 
surface temperature. Effects of lower accommodation co- 
efficients will also be assessed. 


Flow Boundary Conditions 

In simulating highly-rarefied flows, the computational 
flow domain must extend far enough upstream of the body 
to provide ample opportunity for free-stream molecules to 
interact with those molecules that have reflected from the 
body and are diffusing into the flow. Insufficient upstream 
domain size leads to overprediction of aerodynamic heating 
and forces. 2,13 

Taking advantage of body symmetry, the simulated flow- 
field configuration is that of a wind tunnel depicted for pitch 
simulations in Fig. 4 (100 Lx84Hxl8Win cells) with a 
specularly-reflecting symmetry plane and a free-molecular 
outer plane. Particles which strike the outer plane after re- 
flecting off of a body surface are removed from the flow 
field. All other particles striking the outer plane reflect 
specularly back into the flow. This effectively reduces the 
influence that the outer wall has upon the flow near the ve- 
hicle in highly-rarefied or free molecular flows. 

Particles enter the flow domain from the left plane at 
the specified pitch angle, and exit normally from the right 
plane. The top and bottom planes represent free-molecular 
periodic boundary conditions. Again, particles which had 
reflected from the body before striking these planes are 
removed from the simulation. All other particles strik- 
ing one of these planes simply re-enter the flow domain 
from the opposite plane, maintaining their velocity vectors. 


This boundary condition still provides free-stream flow ap- 
proaching the body geometry while reducing the required 
height of the wind tunnel, therefore minimizing the size and 
cost of the simulation. 

SIMULATION RESULTS 

Employing a free-stream flow speed of =8.6 km/s 
and the day-side atmospheric data of Keating 14 , simulations 
were performed for several entry attitudes of pitch, yaw, 
and solar panel canting. Flight conditions are listed in Ta- 
ble 1 for each altitude. Results from the particle simula- 
tion for rarefied flows, employing the finite molecular mean 
free path, were compared to free molecular (collisionless) 
simulation results. These were, in turn, compared to ana- 
lytic free molecular results generated with the FREEMAC 
code 3 which employs ray-tracing algorithms to assess aero- 
dynamic coefficients of spacecraft. FREEMAC results as- 
sumed that all surfaces have a fixed temperature of 300 K. 
The DSMC code predicted surface temperatures near 400 K 
for the solar panels using the radiative equilibrium model 
described above. 

Limited computational resources and the large number 
of cases investigated in the present study restricted the 
number of particles which could be employed in the sim- 
ulation to just four particles per cell in the free-stream. 
However, density gradients in the flow-field significantly 
increased the particle densities near body surfaces such 
that 400,000 particles existed typically in the flow-field 
at steady-state. Employing roughly 6,000 transient steps 
and 6,000 steady-state sampling steps in the simulation, 
the number of statistical samples was sufficiently large to 
yield accurate and meaningful solutions. Run-times aver- 
aged 1. 0-2.0 CPU hours per case on the Cray-YMP su- 
percomputer with a nominal computational performance of 
roughly 1.1 /isec/particle/timestep. The simulation could 
have employed ten times as many particles but would have 
required a ten-fold increase in computational time per case 
which was not warranted in the present study. 

The simulation computed the net force and heat flux 
upon each surface facet of the vehicle. Torques were com- 
puted by the moments of body forces about a reference 
point defined by the intersection of the central body axis 
(x-axis) and the solar panel axis (y-axis). This reference 
point is very close to the spacecraft center of mass. Coef- 
ficients for forces (F) and moments (A/) are defined from 
free-stream properties u ^ and , the simulated frontal 
area A j , and a reference length for moments given by the 
HGA diameter D . 

Pitch and Yaw Aerodynamics 

Moment coefficients for pitch and yaw attitudes (0° - 
90°) are plotted in Figs. 6 and 5 as computed with the 
FREEMAC and DSMC codes. The thermal accommoda- 
tion coefficient was fixed at unity for both codes. For all 
cases, the resulting moments acted to return the vehicle to 
its nominal zero-pitch, zero-yaw attitude. Yaw moments 
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rose smoothly, reaching a peak at roughly 45° before drop- 
ping again. The pitching moment, however, had two local 
maxima (at 20° and 70°) and a local minimum (at 50°). 
This behavior resulted from the solar panels shading the 
HGA partially from the incident flow at those attitudes. 

The DSMC method was used to simulate rarefied flows 
(finite mean free path) at 136 km and 140 km. Repeating 
these simulations under free molecular conditions yielded 
results which were nearly indistinguishable between the 
two altitudes. Note that the DSMC free molecular results 
agree fairly well with those from FREEMAC. Slight differ- 
ences are likely due to different surface temperature mod- 
els employed in each code and the fact that FREEMAC ne- 
glects multiple reflections of molecules with vehicle sur- 
faces. More important, however, is the observation that 
rarefied and free molecular DSMC results do not differ dra- 
matically, indicating that molecular collisions do not alter 
the vehicle aerodynamics significantly. The moments for 
rarefied flow at 136 km are roughly 7-10% below the free 
molecular results, except at high pitch and yaw attitudes. 
Results are similar at 140 km. 

Other aerodynamic properties under pitch and yaw are 
plotted in Figs. 8 and 7. These include force coefficients 
in the direction of the body axis (F x ) t in the perpendicular 
directions ( F y and F z ) t and in the flow direction (drag, C D ). 
The net heat flux q on the windward solar panel is plotted 
relative to the incident free-stream energy flux, Q . For all 
coefficients, the effects of flow rarefaction led to small yet 
noticeable differences from free molecular results. These 
differences are greatest at lower pitch and yaw attitudes, 
where the solar panels and HGA led to greater wetted area 
of the vehicle. At attitudes near 90°, the wetted area was 
significantly smaller such that the flow was effectively more 
rarefied, leading to little distinction between rarefied and 
free molecular results. 

Effect of Body Resolution 

The FREEMAC and DSMC results above employed the 
same coarse body resolution. This was sufficient to demon- 
strate the small differences between rarefied and free molec- 
ular flow models. However, the FREEMAC code may em- 
ploy a fine body resolution which more precisely resem- 
bles the actual vehicle configuration, including small fea- 
tures such as the MGA, LGA, ALTA, and REM. These 
small features would not alter the effects of flow rarefaction 
demonstrated above, but they do impact the magnitudes of 
vehicle aerodynamics. FREEMAC moment and drag co- 
efficients for each body resolution are plotted in Figs. 9 
and 10. Note that the same general pitch and yaw charac- 
teristics result for each body resolution, but that the magni- 
tudes of each differ somewhat Accurate estimation of the 
vehicle aerodynamics would therefore require adjusting the 
fine-resolution FREEMAC results to account for the small 
rarefaction effects demonstrated with the coarse-resolution 
above. 


‘‘Windmill” Roll Aerodynamics 

A creative experiment suggested by Lyons and Hurlbut 15 
for Magellan aerobraking would involve canting the so- 
lar panels in opposite directions to create, a “wind mil l” 
configuration with the spacecraft Measuring the restor- 
ing roll torques induced by the spacecraft control sys- 
tem would provide insight into rarefied hypersonic gas- 
surface interaction. This windmill configuration was sim- 
ulated with the DSMC code in the present work for 
flow in the nominal direction (along the body-axis) at 
138 km altitude. The solar panels were canted at angles 
of 0°, 7.2°, 22.0° ,4 1.4°, 60°, 75.5°, and 90°, where 0° 
is perpendicular to the incident flow. The thermal accom- 
modation coefficient was fixed at unity. Coefficients for 
drag, roll moment, and panel heating are plotted in Fig. 1 1 . 
Drag and heating results exhibit the same differences be- 
tween rarefied and free molecular flows as predicted above 
from the pitch and yaw simulations. However, the roll 
moment coefficients M x in the rarefied simulations were 
nearly twice as large as those in the free molecular flows. 
Although these roll moments are very small, and there- 
fore statistically sensitive, they exhibit this behavior con- 
sistendy for all cant angles. 

This surprising behavior is unlike all other aerodynamic 
observations in this study, leading to concern about solu- 
tion validity. It was suggested that perhaps the DSMC flow 
domain was so small that particles diffusing out of the peri- 
odic simulation boundaries were reflecting back toward the 
body, thus altering the roll characteristics of the vehicle. 
The simulations were repeated three times with consecu- 
tively larger lateral domains (in y and z directions). The dis- 
tance between the vehicle surfaces and the nearest bound- 
ary were increased by factors of three, six, and twenty. 
Nonetheless, the roll moment coefficients for each case dif- 
fered from the original by less than five percent These re- 
sults indicate that, indeed, the roll moment for the windmill 
configuration is small but highly sensitive to molecular col- 
lisions in the rarefied flow-field. 

Effects of Surface Accommodation 

In the DSMC code, a thermal accommodation coeffi- 
cient of, say, a=0.8, means that particles have a 20% chance 
of reflecting specularly from a surface and an 80% chance 
of reflecting diffusely with energies sampled at the surface 
temperature. Windmill roll-moment calculations with pan- 
els canted at 41.4° were repeated using thermal accommo- 
dation coefficients of a = {0.6, 0.7, 0.8, and 0.9}. Re- 
sults for free molecular and rarefied flows are plotted in 
Fig. 12. With increasing thermal accommodation coeffi- 
cient, the roll moment dropped rapidly while panel heating 
increased. More importantly, the distinction between rar- 
efied and free molecular results was significant only for co- 
efficients above a=0.80. This observation further supports 
the assertion that the windmill experiment would provide 
valuable flight data regarding surface accommodation and 
flow-field properties under these conditions. 



5 AIAA Paper No. 93-3676 (Naas & Schmitt) 


CONCLUDING REMARKS 

The aerodynamics of the Magellan spacecraft during 
proposed entry into the atmosphere of Venus may be com- 
puted for free molecular flows with the FREEMAC code. 
However, at nominal altitudes near 140 km, where the 
free-stream Knudsen number drops below Kn=10, the free 
molecular flow assumption could be questionable. This 
study computed the body torques on the spacecraft, when 
entering with various off-design attitudes, using a DSMC 
particle method for rarefied flow conditions. Repeating 
these simulations for collisionless flow provided a direct 
means of assessing the significance of molecular interac- 
tions in the flow, as well as providing a means to validate 
the technique through comparison to FREEMAC results. 
At all pitch and yaw attitudes, molecular collisions tended 
to reduce the aerodynamic forces, moments, and heating 
below free molecular values by 7-10%. All body torques 
acted to restore the vehicle to its nominal zero-pitch, zero- 
yaw attitude. Canting of the solar panels yielded roll mo- 
ments which were highly sensitive to the thermal accom- 
modation coefficient for gas-surface interaction. For coef- 
ficients near unity, the predicted roll moments were notably 
larger for rarefied flows than for free molecular flows. In 
general, however, the DSMC calculations verified that the 
aerodynamics during Magellan aerobraking may be mod- 
eled within 10% error by assuming free molecular flow at 
the nominal entry altitudes. 
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TABLE 1: Simulation Flowfield Specifications 


Altitude 

(km) 

(kg/m 3 ) 

(K) 

m 

(g/mole) 

Woo 

Kn ro (a) 

136 

1.39 x 10~ 8 

217 

40.1 

35.2 

2.819 

138 

9.24 x 10 -9 

222 

39.4 

34.6 

4.222 

140 

6.13 x 10~ 9 

225 

38.9 

34.0 

6.324 


(a) Based on HGA diameter of 370 cm, Knoo = \oo/D. 



Fig. 1 Magellan spacecraft configuration versus simulated ce ^ s * n hzigih (f° r 45 study), and 12 cells in width (from 
geometry shown at nominal attitude (zero pitch and yaw). s y mmetr y plane). 



Location Relative to Solar Panel (Cells) 


Fig. 2 Profiles of number density and molecular mean free 
path along the stagnation streamline ahead of the solar panel 
at its center (nominal attitude and flight conditions). 



Yaw Attitude (Degrees) 

Fig. 5 Coefficient of restoring yaw moment versus en- 
try yaw attitude. Comparison of rarefied DSMC results 
(at 136 km and 140 km) and free molecular DSMC and 
FREEMAC results. 



Fig. 3 Translational temperature profile along the stagna- 
tion streamline ahead of the solar panel at its center (nomi- 
nal attitude and flight conditions). 



Pitch Attitude (Degrees) 

Fig. 6 Coefficient of restoring pitching moment versus entry 
pitch attitude. 
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Yaw Attitude (Degrees) 


Attitude (Degrees) 


Fig. 7 Aerodynamic coefficients versus yaw attitude for Fig. 10 Coefficients for drag versus pitch and yaw attitudes 

drag (C D )> solar panel heating (q/Q)> and forces along body predicted from FREEMAC code using coarse and fine body 

axis ( F x ) and perpendicular axis (F z ). resolutions. 



Pitch Attitude (Degrees) S° ,ar Panel Cam An 9 ,e (Degrees) 

Fig. 8 Aerodynamic coefficients versus pitch attitude for Fig. 11 Coefficients for drag, solar panel heating (q/Q), and 

drag (C D ), solar panel heating ( q /Q) y and forces along body torque (M x , multiplied by 20) versus solar panel cant angle, 

axis ( F x ) and along panel axis ( F y ). Comparison of rarefied results at 138 km to free molecular 

predictions. Surface accommodation fixed at a = 1.0. 



Attitude (Degrees) 

Fig. 9 Coefficients for restoring moments versus entry at- 
titudes in pitch and yaw. Comparison of results predicted 
from FREEMAC code using coarse and fine body resolu- 
tions. 



Solar Panel Thermal Accommodation Coefficient, a 

Fig. 12 Aerodynamic coefficients versus thermal accom- 
modation coefficient a for solar panels canted at 41.4°. 
Comparison of rarefied results at 138 km to free molecu- 
lar predictions. 








